Ab initio Derivation of Low-Energy Model for Alkali-Cluster-Loaded Sodalites 
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An effective low-energy model describing magnetic properties of alkali-cluster-loaded sodalites is derived by 
ab initio downfolding. We start with constructing an extended Hubbard model for maximally localized Wannier 
functions. Ab initio screened Coulomb and exchange interactions are calculated by constrained random phase 
approximation. We find that the system resides in the strong coupling regime and thus the Heisenberg model 
is derived as a low-energy model of the extended Hubbard model. We obtain antiferromagnetic couplings ~ 
O(10 K), being consistent with the experimental temperature dependence of the spin susceptibility. Importance 
of considering the screening effect in the derivation of the extended Hubbard model is discussed. 
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I. INTRODUCTION 

Zeolites constitute a huge family of nonporous crystalline 
aluminosilicates which have a wide variety of intriguing 
properties. 1 Because of their capability of hosting various 
ions, atoms, molecules and clusters in their subnanometric 
pores with rich possibilities of crystal structure, they have 
versatile functionalities such as high catalytic activities, sorp- 
tion characteristics, ion-exchange abilities. Numerous zeolites 
with various compositions and framework topologies have 
been synthesized and exploited in many applications. Be- 
sides such diverse fascinating aspects, it is of great interest 
to focus on electron correlations in this system. A variety 
of zeolites have been viewed as three-dimensional correlated 
i-electron systems providing a non-trivial play ground for a 
systematic control of many-body correlation effects.— For ex- 
ample, although all the ingredients of zeolites are nonmag- 
netic elements, some of zeolites exhibit an intriguing mag- 
netism for certain conditions; zeolites LTA and LSX with 
potassium clusters have ferromagnetic ground states depend- 
ing on the number of potassium atoms per cage, 3,4,5 while 
sodalites loading various alkali-metal clusters exhibit robust 
antiferromagnetisir h 6,7 ! 8 ' 9 

When we study such characteristic many-body effects in 
zeolitic materials, it is definitely impractical to calculate ev- 
erything from first principles. The unit cell is extremely huge 
and contains many atoms, so that formidable numerical cost 
would be required. On the other hand, recent conventional ab 
initio studies have clarified that some zeolites have quite sim- 



ple low-energy electronic structures 
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for example, in 



the sodalite system, the aluminosilicate cage forms a wide gap 
more than 5 eV around the Fermi level and, in this energy 
gap, electronic states due to guest alkali clusters make narrow 
bands with the width ~1 eV. Its band dispersion is well rep- 
resented by simple tight-binding models, which indicates that 
the so-called "superatom" picture 2 or the "particle-in-a-box" 
model 10 correctly captures essential aspects of the low-energy 
physics of the sodalite systems. With this situation, rather 
than the full ab initio approach, the three-stage approach is 
expected to work more successfully. Here, in the first stage, 
we perform standard ab initio density-functional calculations 
and, in the next step, we perform downfolding procedures, 



that is, construction of an effective low-energy model. Finally, 
we solve the resulting model by high accurate and reliable 
solvers. The so-called "LDA+DMFT" method 14 combining 
local density approximation (LDA) and dynamical mean field 
theory (DMFT) is a typical example of this approach. Re- 
cently, the three-stage approach has been extensively applied 
to various correlated electron systems. Especially, it has been 
demonstrated that the scheme really works with high accuracy 
for various transition-metal oxides J^ 5 - 

Recently, as a reliable tool for evaluating the values of inter- 
action parameters in the downfolding step, a constrained ran- 
dom phase approximation (cRPA) method is formulated 
Compared to the standard method based on a constrained 
LDA technique, 18 the cRPA has several advantages; one can 
precisely exclude screening processes between the basis states 
of the effective model, which should be considered in the last 
stage solving the effective model. In addition, we can ob- 
tain onsite and off site interactions at one time. While the 
cRPA method has been employed in many studieS j 15 i 19 ' 20 ' 21 
applications to zeolitic materials have yet to be done. Indeed, 
it is a highly non-trivial issue to determine the values of in- 
teraction parameters of zeolites; the bases of the low-energy 
model of these materials are no longer localized at some spe- 
cific atoms and are extended spatially beyond several guest 
atoms in the cage. So, we have to evaluate the value of in- 
teraction parameters not for atomic orbitals but for molecular 
orbitals. In order to construct automatically such basis func- 
tions with non-trivial spatial spread, it is convenient to exploit 
maximally localized Wannier orbitals (MLWOs). 22 Recently, 
MLWO is combined with cRPA calculations to estimate the 
onsite Hubbard U as well as offsite interaction parameters in 
the low-energy models of various systems. 2 - ^ 

It is of great interest to apply this state-of-the-art downfold- 
ing technique based on the combination of cRPA with MLWO 
to zeolitic materials and examine how it works. As a bench- 
mark for this purpose, we consider sodalites which are clas- 
sified as the simplest zeolite. The framework of this mate- 
rial is described as a body-centered-cubic (bcc) array of f3 
cages [(Si02)3(A10^)3, the smallest unit of the aluminosil- 
icate cage] and each cage accommodates ionic alkali clus- 
ters A^ + to neutralize the negative charge of the framework. 
Experimentally, it has been well established that the system 
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can be viewed as a crystal of F centers sitting on the center 
of the tetragonal cluster A| + . Especially, a magnetic prop- 
erty of a sodium electrosodalite (or black sodalite) is quite 
well understood in terms of the 5=1/2 Heisenberg model on 
the bcc latticei 6 ' 7 ' 8 ' 9 Since the temperature dependence of the 
magnetic susceptibility of the Heisenberg model can be cal- 
culated by the high-temperature-expansion scheme with high 
accuracy, we can obtain the precise values of exchange cou- 
pling by parameter fitting to the experimental data. Thus, the 
sodalites are the best systems to examine the reliability of de- 
rived parameters with the ab initio downfolding, in that we 
can compare unambiguously the theoretical exchange values 
and the experimental ones. 

The purpose of the present study is to examine how accu- 
rately we can construct a low-energy model of the sodalite 
system by the ab initio downfolding technique. While one 
can exploit direct ab initio calculations based on local spin 
density approximation to evaluate the exchange coupling, 23 
the present study focuses on an alternative approach which 
is feasible to not only localized spin systems but also more 
general cases. First, we construct a single-band extended 
Hubbard model based on the ab initio downfolding scheme 
and then derive an effective Heisenberg model by the second- 
order perturbation. A similar strategy was taken in the previ- 
ous study 10 but, there, the "kinetic-exchange" term 24 only was 
evaluated and the "direct-exchange" term was completely ne- 
glected. As shown below, in the sodalite, the direct exchange 
has the same energy scale as the kinetic exchange and thus the 
two exchange couplings compete with each other. In addition, 
in the past parameter estimations, the screening effect was 
completely neglected. We will show the importance of tak- 
ing the screening effect into account in the parameter deriva- 
tion; if we neglect the screening effect, the kinetic-exchange 
value is smaller than the direct-exchange value and the net 
exchange becomes ferromagnetic. When the screening effect 
is switched on, the kinetic exchange reverses the direct ex- 
change, thus resulting in antiferromagnetic interactions being 
consistent with the experiments. 

This paper is organized as follows. In Sec. II, we de- 
scribe our basic strategy for deriving the effective Heisenberg 
model from first principles. Section III is devoted to the accu- 
rate estimation of the experimental exchanges using the high- 
temperature expansion to the Heisenberg model. Following 
by recent measurement of the magnetic susceptibility for the 
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and 

Jij = {4>i4>j\W\<f)j(f>i) 



where W(r,r') is a screened Coulomb interaction. Vij at i = 
j corresponds to onsite Hubbard parameter U. 



sodalites, 9 we give the exchange parameters of the sodium 
electrosodalite and the potassium electrosodalite and discuss 
the differences between the two. In Sec. IV, ab initio compu- 
tational results are presented and compared with the experi- 
mental results. The concluding remarks are given in Sec. V. 



II. AB INITIO CONSTRUCTION OF EFFECTIVE 
HAMILTONIANS 

We consider ab initio derivations of the effective Heisen- 
berg model describing "low-energy" electronic structures. 
Conventionally, the derivation is based on the second-order 
perturbation to the single-band extended Hubbard Hamil- 
tonian consisting of the transfer part Ht, the Coulomb- 
interaction part TLv, and the exchange-interaction part Hj as 

H = H t + H v + Hj (1) 

with 

a ij 



op ij 



crp ij 

where a\ a (a,i a ) is a creation (annihilation) operator of an elec- 
tron with spin a in the Wannier orbital localized in the ith so- 
dalite cage. The parameters in Eq. (f2]i contain an onsite 
energy (i = j) and hopping integrals (i ^ j), written by 

Uj = (&|Wo|&) (5) 

with \<j>i) = aAO) and Ho being the one-body part of H. The 
Vij and Jij values in Eqs. (f3]l and (01 are screened Coulomb 
and exchange integrals in the Wannier orbital, respectively, 
expressed as 



(6) 



(7) 
I 

Now, we consider a situation with the half-filling and 
atomic-limit condition, where the parameters satisfy the fol- 




drdr'4>* (r)4>i (r)W(r, r')# (r% (r') 
drdr'cj)* (r)<pj (r)W{r, r>* (r')& (r'), 
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lowing inequality 



U-Vtj > \tij\ > o. 



(8) 



In this situation, with the second-order perturbation, the effec- 
tive Hamiltonian which describes the fine "low-energy" spec- 
trum associated with the spin structure is given as the follow- 
ing Heisenberg model^ 5 . 



Ti-cS — 2 SfijSi ■ S 

i>3 



31 



(9) 



where the local spin operator Si is conventionally represented 
in term of the creation and annihilation operators as Sf — 

\{a\] a H + a li a i^ S? = Yi( a lt a H ~ a ll a *T)> and S i = 
5(0]+ ciif — ot,aij_). The effective exchange coupling in Eq. (0 
is written as 



jij — Kij J, 



with 



2\tj£ 
U - Va 



(10) 



(11) 



The Kij is a "kinetic-exchange" term 24 which stabilizes the 
antiferromagnetic coupling between the local spins, while the 
second term in Eq. ( TTOb is a "direct-exchange" term favor- 
ing the ferromagnetic coupling. The competition between the 
two-type exchange terms determines the net magnetic feature 
of the system (i.e., whether the system prefers the antiferro- 
magnetic or ferromagnetic state). 

The calculation of the effective exchanges Jij is basically 
straightforward after parameterizations of t y -, Vij, and Jy 
but a careful treatment is needed for the calculation of the 
screened Coulomb interaction of W(r, r'). The screened in- 
teraction considered in the extended Hubbard model should 
not include screening formed in a target band of the model. 
This screening should be considered at the step of solving 
the effective model and, at the downfolding stage, we must 
exclude the target-band screening effects to avoid the double 
counting of this screening. In the random phase approxima- 
tion (RPA), this constraint is easily imposed, 16,17 because the 
RPA polarization function is given as the sum of the band pairs 
associated with individual transitions; we first calculate the 
polarization function with excluding the transitions in the tar- 
get band and then evaluate the screened interaction W(r, r') 
with using this polarization function. Finally, we compute the 
Vij and parameters as the Wannier matrix elements of the 
W interaction. 

There are two other choices on the treatment of the 
Coulomb interaction. The first is the use of "bare" Coulomb 
interaction v(r, r ') = ]r T ~r 7 J mstea( l °f W(r, r'). The result- 
ing Vij and Jij parameters have no screening effect and will 
give larger values than the constrained-RPA values discussed 
above. The kinetic-exchange parameter Kij becomes small 
because of the increase of U—Vij in Eq. (fTTT i. We note that this 
choice has been widely used in the literature so fa r 10 ! 11 ' 12 ! 13 
but there is no justification. Another choice is the use of the 



"fully" screened Coulomb interaction, where we calculate the 
RPA polarization function with no constraint on the transi- 
tions. The result includes the target-band screening effect and 
therefore the calculated Vij and Jy values will be largely re- 
duced compared to the constrained-RPA values. We compute 
the interaction parameters Vij and Jij with the "bare", "con- 
strained RPA", and "full RPA" interactions and discuss the 
importance of the screening effect on the derivation for ex- 
change values of the Heisenberg model. 



HI. ESTIMATION OF THE EXCHANGE COUPLINGS 
FROM THE EXPERIMENT 

Before presenting ab initio computational results, we con- 
sider experimental values of exchange couplings, which are 
estimated from the data for the temperature dependence of 
the magnetic susceptibility. In the sodium electrosodalite, 
the measured Weiss temperature 0w is —170 K, while the 
Neel temperature Tn is 50 K. The negative Weiss tempera- 
ture and the existence of the antiferromagnetic transition in- 
dicate the antiferromagnetic interaction between neighboring 
spins, while inequality |0w|3>In implies that there is strong 
frustration in the system or equivalently the presence of next- 
nearest-neighbor exchange couplings. In fact, magnetic prop- 
erties of the sodium electrosodalite have been discussed with 
the Heisenberg model up to the next nearest neighbors. Re- 
cently, the magnetic measurement has been performed for the 
potassium electrosodalite and the 0w and Tn temperatures 
are observed as —330 K and 80 K, respectively^ Here, we 
determine the exchange parameters in the Heisenberg model 
so that the calculated model Weiss and Neel temperatures re- 
produce the experimental ones. Accuracy of 9w and Tn ob- 
tained from solving the model critically affects the quality of 
the exchange couplings. In this work, we calculate the high- 
temperature series of the spin susceptibility up to tenth order 
in inverse temperature 35 using the finite cluster method^ 

The explicit form of the Heisenberg model on the bcc lattice 
up to the next nearest neighbors is given as 



H. 



err 



2J X Si 

(ij) 



Sj 



2j 2 s * 

(ay 



(12) 



where the first summation (ij) is taken over the bonds be- 
tween nearest neighbors and the second summation (ij) over 
the bonds between next-nearest neighbors. J\ and J<z rep- 
resent the exchange couplings for the nearest neighbors and 
next-nearest neighbors, respectively. Note that the suffices 
"1" and "2" attached to J hereafter specify the bond between 
the nearest neighbors and the bond between the next-nearest 
neighbors, respectively. 

The spin susceptibility for a general wavevector, x(g), can 
be expressed as^ 5 - 



1 ft 
X{q) =Nj 



Sf 



S*j)e 



(13) 
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Here, j3 is the inverse temperature and 
resents the thermodynamic average; i.e., 



rep- 



Tr(- 



-PHefi 



Uniform and staggered 



spin susceptibilities are given as x = x(0) and x(Q) with 
Q = (it, it, it), respectively. The \ and x{Q) U P to the m " st 
order in (3 are given as 



A X T= 1-(3(4J 1+ 3J 2 ) + 0(/3 2 ), 
4 X (Q)T = 1 + (3{AJ X - 3J 2 ) + 0(f3 2 ). 



(14) 
(15) 



The first-order coefficients — (4Ji+3j72) and4j7i — 3J72 above 
correspond to the high-temperature-limit Weiss temperature 37 
and the mean-field Neel temperature, 38 respectively. It should 
be noted here that the temperature range for experimental 0w 
is far from the high-temperature limit, so that the Weiss tem- 
perature given above is not a good estimate for the experi- 
mental value. Furthermore, the mean-field Neel temperature 
is seriously overestimated because the quantum fluctuation is 
neglected. 

To go beyond the first-order analysis and obtain the pre- 
cise temperature dependence of the magnetic susceptibility, 
we consider the higher-order expansion series of \ and x(Q)- 
In order to extrapolate the series down to low temperatures, 
we use the Pade approximation, in which a series of x and 
x{Q) is approximated as tJ-^jtj, where Pl{x) and Qm(x) 
are the i-order and Af-order polynomials, respectively. We 
call it the [L, Al] Pade approximation. Figure Q] illustrates an 
example of the extrapolation of x 1 at Ji = 0. We find that 
various Pade approximations show good convergence down to 
T/Jx~3. To estimate the Weiss temperature, 0w, we fit x as 
X = (T — 0w) _1 as shown in the dotted line in Fig.Q] The 
fitting temperature range is 5<T/j7i<10. (As shown below, 
this range roughly corresponds to the experimental tempera- 
ture range.) The 9w estimated as — 6.OJ1 indicated by the 
arrow in the figure is appreciably smaller than the first-order 
value — AJ\. The treatment can straightforwardly be applied 
to the case of J2^0. In Fig. [2] the calculated Weiss tempera- 
ture (dots) is shown as a function of Jij J\. 

The estimation for Tn proceeds as follows: x(Q) is ex- 
pected to have a pole at the finite (3 and behaves as 



X (Q)cx(/?-/3 N )-\ 



(16) 



where /?n is the inverse of the antiferromagnetic transition 
temperature (Tn = 1/ /?n) and the 7 is the critical exponent. 
In the three-dimensional Heisenberg model, 7 is known to 
be ~1 .39 (Ref. [39h. By taking the logarithmic derivative of 
Eq. (TToT l. we obtain 3 ^ 



rfiogx( Q) 

d(5 



x 



"7 



(17) 



Since the Pade approximation can describe simple poles ex- 
actly, approximations to the logarithmic derivative should 
converge much faster. In addition, we can evaluate 7 from a 
residue of the pole as well as the location of the pole giving the 
critical temperature. At J2 = 0, various Pade approximations 
show good convergence and give = 2.76 and 7 = 1.39. 
For finite J2, however, the frustration lowers the Neel temper- 
ature and it becomes difficult to estimate Tn and 7 accurately. 
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FIG. 1: (Color online) Inverse of the uniform spin susceptibility at 
J 2 = 0. [3,3], [4,4], [5,5], [4,5], and [5,4] Pade approximations 
are shown. We also plot the Curie-Weiss fitting as the dotted line. 
An arrow indicates the extrapolated Weiss temperature estimated as 
—6.0^71. Notice that this value is largely deviated from the high- 
temperature-limit value —4J\. 
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FIG. 2: Weiss temperature, 0w, and Neel temperature, Tn, as a 
function of Jil J\. Vertical errorbars for Tn represent the scatter- 
ing of the various Pade approximation. For sodium electrosodalite, 
the experimental 9w (—170 K) and Tn (50 K) are well reproduced, 
when J\ and J2I J\ are set to 26 K and 0.31 (vertical dotted line), 
respectively. In the case of potassium electrosodalite, such J\ and 
Jil J\ were found to be 48 K and 0.42 (vertical dashed line), re- 
spectively. 



Thus, to improve the convergence, we use the Pade approxi- 
mation of x 1 ^ 7 with 7 kept at 1.39. This assumption works 
well even for the finite J2 and we obtain Tn as a function 
of J2IJ1 (crosses of Fig. |2). The errorbar comes from the 
scattering of the various Pade approximations. By using these 
data and referring the experimental 0\y and Tn temperatures, 
we reasonably estimate the exchange couplings of the sodalite 
system. The resulting values are J\ = 26 K and J 2 = 8 K for 
the sodium electrosodalite and J\ = 48 K and J2 = 20 K for 
the potassium electrosodalite. 
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IV. RESULTS AND DISCUSSIONS 

Our ab initio calculations were performed with Tokyo Ab 
initio Program Packaged With this program, electronic- 
structure calculations with the generalized-gradient- 
approximation (GGA) exchange-correlation functional 29 
were performed using a plane-wave basis set and the 
Troullier-Martins norm-conserving pseudopotentials 30 in the 
Kleinman-Bylander representation. 31 ' 32 The energy cutoff in 
the band calculation was set to 49 Ry and a 5x5x5 /c -point 
sampling was employed. The experimental crystal-structure 
data were taken from Ref. for sodium electrosodalite 
and Ref. [l~2] for potassium electrosodalite. The calculations 
for the screened interactions are followed by Ref. |2l|. The 
polarization function was expanded in plane waves with 
an energy cutoff of 5 Ry and the total number of bands 
considered in the polarization calculation was set to 200. The 
Brillouin-zone integral on wavevector was evaluated by the 
generalized tetrahedron method. 33 The additional terms in the 
long-wavelength polarization function due to nonlocal terms 
in the pseudopotentials were explicitly considered following 
Ref. [34]. A problem due to the singularity in the Coulomb 
interaction, in the evaluation of the Wannier matrix elements 
Vij and Jij, was treated in the manner described in Ref. [34]. 

We show in Fig.[3]ab initio GGA band structures (red solid 
lines) of (a) sodium electrosodalite and (b) potassium elec- 
trosodalite. We see an isolated band near the Fermi level (en- 
ergy zero). This band is due to confined electrons in the so- 
dalite cage and we employ this band as the target band of the 
extended Hubbard model. The entangled band structures be- 
low —4 eV and above +1 eV are associated with electronic 
states of the framework of the sodalite. The overall band 
structure of the sodium electrosodalite is similar to that of the 
potassium electrosodalite. A notable difference is that the tar- 
get bandwidth of the sodium electrosodalite is 0.86 eV, while 
that of the potassium electrosodalite is 1.01 eV, which makes 
differences in the values of transfer integrals of the two mate- 
rials (see below). 

Figure [4] visualizes our calculated maximally localized 
Wannier orbitals for the target band of (a) sodium electroso- 
dalite and (b) potassium electrosodalite. We can see that the 
resulting Wannier orbital is confined in the cage and has an 
s symmetry around the cage center. The calculated spatial 
spread of the Wannier orbitals are 2.66 A for the sodium elec- 
trosodalite and 2.91 A for the potassium electrosodalite and 
these values are smaller than the diameter of the cage (7.6 
A for the sodium electrosodalite and 8.0 A for the potassium 
electrosodalite). 

We next calculate transfer integrals in Eq. <(5j as matrix 
elements of the Kohn-Sham Hamiltonian in the Wannier or- 
bital. The nearest-neighbor transfer t-y and the next-nearest- 
neighbor transfer t2 are found to be —57.3 meV and —32.1 
meV for the sodium electrosodalite. The results for the potas- 
sium electrosodalite are —68.0 meV and —31.1 meV. It was 
found to be negligibly small for other transfers beyond the 
third neighbors; their magnitudes are less than a few meV. The 
band dispersion calculated with t\ and i 2 is shown as blue dots 
in Fig. [3] We can see that the original band structure is quite 




FIG. 3: (Color online) Calculated ab initio band structures (red solid 
lines) of (a) sodium electrosodalite and (b) potassium electrosodalite. 
The blue dotted dispersions are obtained by the t\-t2 model, where 
ti and t2 are nearest and next-nearest transfers, respectively. For the 
values, see the text. The zero of energy is the Fermi level. 



well reproduced with the two-parameter model. We note that 
the Kohn-Sham Hamiltonian Hks is different from the exact 
one-body Hamiltonian Ho in Eq. ©. The difference between 
the two requires involved discussions about the "downfolding 
self energy"-l&i2 so, in the present study, for the simplicity, 
we employed the Hks instead of the Hq. 

Figure[5]plots the Wannier matrix elements of the screened 
Coulomb interaction Vy (green dots) calculated with con- 
strained RPA, as a function of the distance between the cen- 
ters of the MLWOs; r = \{4>i\r\(j>i) - (<t>j\r\<t>j)\. The pan- 
els (a) and (b) show the results of the sodium electrosodalite 
and of the potassium electrosodalite, respectively. The Vij 
decays as an isotropic function of l/(er) (dotted line) where 
e is a macroscopic dielectric constant calculated with cRPA. 
The value of e is 3.2 for the sodium electrosodalite and 3.0 for 
the potassium electrosodalite. For comparison, we also plot 
bare Coulomb interactions (red dots), which should decay as 
1/r (solid line) beyond the nearest-neighbor distance (>7 A). 
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FIG. 4: (Color online) Calculated maximally localized Wannier 
functions of (a) sodium electrosodalite and (b) potassium electroso- 
dalite. The amplitudes of the contour surface are +1.5/\/v (blue) and 
— 1.5/i/v (red), where v is the volume of the primitive cell. Si, O, Al, 
and Na or K nuclei are illustrated by blue, silver, yellow, and green 
spheres, respectively. 



We see that the bare Coulomb interaction is reduced in less 
than half by considering the screening effect with cRPA. On 
the top of this, the full RPA screened Coulomb interactions 
are shown as blue dots, which are nearly zero, except for the 
onsite value at r = 0. The exchange interactions of Jy were 
found to decay very quickly; the magnitude is nearly zero, ex- 
cept for the nearest and next-nearest values. This quick decay 
was the same for the three cases of the bare, cRPA, and full 
RPA. 

We summarize in TABLE U the principal parameters in the 
resulting Hubbard model of Eqs. (d}-©; the interactions up 
to the next nearest neighbors. The table compares U, V±, V2, 
J 1, and J r 2 calculated with the three-type interactions of the 
bare, cRPA, and full RPA. We see that the calculated values 
become small in order of increasing the screening (bare — > 
cRPA — > RPA). From cRPA to full RPA, the diagonal parts 
of the Coulomb interaction, U, V\, and V2, are significantly 
reduced by more than one order of magnitude. This is simply 
because a metallic screening process is switched on at the full 
RPA. In contrast, the off-diagonal parts of the Coulomb inter- 
action, Ji and J2, are not so screened. 20 The kinetic-exchange 
values K\ and K2 obtained via Eq. (fTTT i. in the bottom two in 
the table, exhibit increasing tendency with bare — > cRPA — > 




5 10 15 20 25 30 
r (Angstrom) 




5 10 15 20 25 30 
t (Angstrom) 



FIG. 5: (Color online) Calculated screened Coulomb interactions 
of (a) sodium electrosodalite and (b) potassium electrosodalite as a 
function of the distance between the centers of maximally localized 
Wannier orbitals displayed in Fig. [4] The red, green, and blue dots 
represent the result with the bare, constrained random phase approx- 
imation (cRPA), and full-RPA interactions, respectively. The solid 
and dotted curves denotes 1/r and l/(er), respectively, where e is a 
macroscopic dielectric constant calculated with the cRPA. 



RPA. 

There are discernible differences between the sodium elec- 
trosodalite and the potassium electrosodalite; for example, for 
cRPA, U of the sodium electrosodalite (2.71 eV) is somewhat 
larger than that of the potassium electrosodalite (2.47 eV) and 
Ji of the sodium case (27.0 K) is nearly half of the potassium 
case (44.5 K). These results can be consistently understood 
in terms of the smaller spatial spread of maximally localized 
Wannier orbitals of the sodium case (2.66 A) than that of the 
potassium case (2.91 A). 

The accuracy of the second-order perturbation in deriving 
the exchange parameters in the Heisenberg model is checked 
by an estimate of (U — Vi)/tx. The values are large enough 
for the both sodalites (in cRPA, 36.7 for the sodium electroso- 
dalite and 27.1 for the potassium electrosodalite), indicating 
that the system is close to the atomic limit enough and the per- 
turbation treatment is reasonably justified. All odd-order con- 
tributions with respect to transfer TL t of Eq. to the kinetic 
exchanges vanish, independent of the lattice. 40 The first cor- 
rection to the second-order perturbation arises from the fourth 
order, which is negligibly small. 

We show in TABLE [TT] the theoretical Heisenberg ex- 
changes Ji and J2 [Eq. (fTOb l obtained with the Hubbard- 
model parameters in TABLE Q] and compare those with the 
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TABLE I: List of interaction parameters in the single-band ex- 
tended Hubbard model in Eq. (Q3, together with kinetic exchanges in 
Eq. Jilt . The parameters obtained with the bare, constrained random 
phase approximation (cRPA), and full RPA are compared. Suffices 
"1" and "2" attached to V, J, and K specify the nearest neighbors 
and the next-nearest neighbors, respectively. Units are eV for U, Vi, 
and Vi and K for Ji, J2, Ki, and if 2. 





sodium electrosodalite 


potassium electrosodalite 




bare 


cRPA 


RPA 


bare 


cRPA 


RPA 


u 


5.79 


2.71 


0.13 


5.34 


2.47 


0.17 


Vi 


1.77 


0.61 


0.01 


1.70 


0.63 


0.01 


v 2 


1.54 


0.54 


0.00 


1.47 


0.54 


0.00 


Jl 


56.9 


27.0 


22.0 


97.2 


44.5 


39.0 


J 2 


22.4 


10.6 


8.6 


20.9 


9.9 


8.1 


K x 


18.9 


36.3 


596.5 


29.4 


58.2 


631.8 


K 2 


5.6 


11.0 


182.5 


5.8 


11.7 


129.2 



experimental values derived in Section III. We see that the 
exchange couplings qualitatively change by considering the 
screening effect; the sign of the couplings changes from nega- 
tive (ferromagnetic interaction) to positive (antiferromagnetic 
interaction) between the bare and cRPA. The values further 
enhance as proceeding from cRPA to full RPA but the latter 
gives a clear overestimate due to the large size of kinetic ex- 
changes (see TABLE|I]i. For the agreement between the theory 
and experiment, the cRPA is clearly the best among the three 
cases of the bare, cRPA, and full RPA. 

TABLE II: List of parameters of the Heisenberg model in Eq. l[9}, 
where J\ and J2 are the nearest and next-nearest exchange cou- 
plings. The theoretical values with the bare, constrained random 
phase approximation (cRPA), and full RPA are compared with the 
experimental results obtained in Sec. III. The unit is K 





sodium electrosodalite 


potassium electrosodalite 




bare cRPA RPA Expt. 


bare cRPA RPA Expt. 


Jl 


-37.9 9.3 574.5 26 


-67.8 13.8 592.9 48 


Ji 


-16.8 0.4 173.9 8 


-15.1 1.8 121.1 20 



However, the calculated values of J with cRPA are still 
quantitatively underestimated from the experiment. This may 
be partially attributed to errors in the derived Hubbard-model 
parameters. A possible error is underestimation of the transfer 
parameters calculated as matrix elements of the Kohn-Sham 
Hamiltonian Hks- The Hks already includes the self-energy 
effect due to electron-electron interactions in the target bands 
of the Hubbard model as the exchange-correlation potential. 
As mentioned above, in the downfolding schemer&ii this self 
energy must be excluded in the stage of the derivation of the 
Hubbard model. If we use the exact Ho not including the 
target-band self energy, the magnitudes of the evaluated trans- 
fers will become larger quantitatively. We found that an ar- 
tificial enhancement of the transfers by 20 % leads to a sat- 



isfactory improvement to the underestimation of J observed 
above. (For example, for the sodium electrosodalite, with this 
modification, the J\ and J2 values change from 9.3 K and 0.4 
K to 25.3 K and 5.2 K, respectively.) 

Another possibility of the error might arise in the inter- 
action parameters evaluated by the constrained RPA. The 
RPA leaves out the vertex correction in the polarization func- 
tion. There are some studies in which the vertex correc- 
tion is treated within local density approximation in density- 
functional frameworki^iiL By considering this effect, the 
screening becomes larger. If we calculate the screened 
Coulomb interaction W(r, r') with using this LDA dielec- 
tric function instead of the RPA one, we will obtain smaller 
values of the interaction parameters. We found that the use of 
an artificially smaller U by 25 % in a J estimation leads to 
an improvement; for sodium electrosodalite, we obtained J\ 
= 26.6 K and J% = 5.4 K. The quantitative discussions about 
the beyond RPA treatment are, however, not simple and need 
to be given more carefully in future studies. 

Finally, we consider an effect of electron-lattice coupling 
on the results. If an electron occupying a superatom local- 
ized s orbital (see Fig. 2]l is transferred to the next site, one 
may then expect relaxation with an orbital expansion, lead- 
ing to a reduction of the onsite Coulomb repulsion; with 
this expansion of the localized orbital, the excitation energy 
to the doubly-occupied state is reduced from U — Vij to 
U — — AS, where AS* is a stabilization energy due to the 
orbital expansion induced by a lattice deformation of tetrahe- 
dral cluster Na4 confined in a (3 cage. Its energy scale can 
be the order ^0. 1 eVr^ and thus taking into consideration of 
this effect is expected to give a substantial improvement. The 
quantitative estimation of AS" from first principles is, how- 
ever, not so easy, which would require to solve technical issues 
including ab initio calculation for electron-lattice coupling. 



V. CONCLUSIONS 

To conclude, we have presented effective Heisenberg 
models describing the magnetic properties of alkali-cluster- 
loaded-sodalite systems. The derivation of the exchange cou- 
plings is based on the second-order perturbation to a single- 
band extended Hubbard model parameterized by ab initio 
density-functional and constrained RPA calculations. Main 
results in the present study are that (i) the direct-exchange 
couplings, dropped in the past studies, were estimated from 
first principles and were found to be the same energy scale as 
the kinetic exchanges and (ii) importance of considering the 
screening effect in the parameter derivation was found out; 
when the screening is properly considered, the net exchange 
couplings Ji and J2 become antiferromagnetic and the re- 
sulting exchange values are in a reasonable agreement with 
the experimental values on the order of ten K. 

In this work, we have considered a single-band system; the 
low-energy electronic structures of the sodalite systems were 
captured in view of a superatom-s-electron picture. It is in- 
teresting to apply the strategy presented here to other zeolites; 
for example, zeolites LTA and LSX described by multi-band 
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systems. The cage size of these materials is bigger than that 
of the sodalite and many alkali atoms more than four can eas- 
ily be doped. As a result, these materials will form partially- 
filled p-band structures of the superatoms. The minimal model 
of these systems are clearly the multi-band model, thus lead- 
ing to a new intriguing magnetic property due to the Hund's 
rule coupling and/or its competition with the kinetic and di- 
rect exchanges. In fact, the temperature-dependence data of 
the magnetic susceptibility of the LTA zeolite strongly suggest 
the possibility of the highly nontrivial ferromagnetic ground 
state.- (The antiferromagnetic behavior suddenly changes to 
the ferromagnetic behavior at 50 K.) There are active debates 
on this mechanism and ab initio calculations aiming at the 
construction of the effective Hamiltonians describing the low- 
energy physics of these systems will helpfully be contributed, 
which remains as future study. 
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